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theoretical  [4,12]  and  computational  advances  relating  to  the  gravity 
model  of  spatial  interaction.  This  paper  details  our  computational 
experience  with  travel  estimation  within  Corpus  Christi,  Texas,  using 
a  new  convex  programming  approach  of  Charnes,  Raike  and  Bettinger  [4). 
It  is  found  that  available  estimation  techniques  necessarily  result  in 
non- integer  solutions.  A  mathematical  device  is  therefore  introduced 
which  assures  integer-values  estimates.  Implications  and  extensions  of 
the  work  are  indicated. 
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INTRODUCTION 

As  reported  elsewhere  [3,  6,  11]  a  large  scale  modelling  effort 
is  underway  for  the  purpose  of  developing  land  use  and  environmental 
guidelines  in  the  Texas  Coastal  Zone.  The  land  use  portion  of  this  pro¬ 
ject  requires  the  spatial  allocation  of  commercial,  retail  and  residential 
elements  of  the  Corpus  Christi  region  under  a  series  of  policy  constraints. 
As  an  input  to  this  analysis  it  has  been  important  to  determine  the  character 
of  spatial  interaction  for  a  series  of  activities  not  the  least  of  which  was 
the  intra- urban  traffic  pattern.  A  set  of  experiment?  have  been  conducted 
on  various  interaction  estimation  procedures  and  the  results  of  one  of 
these  experiments  is  reported  below.  The  importance  of  this  experiment 
lies  in  the  demonstrated  applicability  of  an  extremal  method  of  statistical 
estimation  based  on  the  statistical  information-theoretic  approach  outlined 
by  Charnes,  Raike  and  Bettinger  (4).  The  demonstration  rests  on  a  new 
mathematical  device,  herein  introduced  for  the  first  time.  This  device 
incorporates  our  desire  for  integer -valued  estimates  of  the  traffic  flows, 
and  in  fact  guarantees  integer  solutions.  Past  literature  and  computation 
has  failed  to  come  to  grips  explicitly  with  this  problem.  We  believe  this 
to  be  the  first  full  scale  application  of  the  approach  and  feel  that  some  of 
the  issues  faced  may  be  of  value  to  other  spatial  analysts  concerned  with 


interaction  estimation  problems. 


I.  THE  ESTIMATION  PROBLEM 

The  Texas  Highway  Department  fl3]  has  compiled  and  up-dated 
origin-destination  (O-D)  data  for  the  fifty-six  traffic  districts  within 
Corpus  Christi,  both  for  vehicle  trips  and  for  passenger  trips,  with  each 
trip  type  disaggregated  by  vehicle  type.  Summing  these  data  yielded  the 
aggregate  trip  type  "TOTAL  INTERNAL  PERSON  TRIPS-ALL  MODES  OF 
TRAVEL.  "  We  considered  our  first  estimation  task  to  be  the  district- 
to-district  trip  distribution  of  this  trip  type. 

We  wished  to  estimate  a  set  {t.  :  i,  j  =  l,  .  .  .  ,  56 1,  where  t^  is 

the  number  of  persons  i raveling  from  district  i  to  district  j  in  a  day's 

time.  The  t- ■  were  to  satisfy  a  "gravity  equation"  of  the  form 
J 

ri°) 

tij  -  - —  (1) 

where  K.^  is  the  "generalized  distance"  between  districts  i  and  j,  t.  and 
o  (i,  i=l,  ....  56)  are  "gravity  potentials",  and  g  is  a  function  whose  form 
is  determined  by  the  estimation  procedure.  The  O-D  data  and  the  are 
known  constants,  and  are  inputs  to  the  estimation  procedures.  In  the  present 
effort,  the  Ky's  are  a  function  of  the  road  distance  between  districts  and 
of  the  speed  limits  on  these  roads;  they  were  approximated  by  inspection 
of  maps  of  the  city. 

The  trip  distr  ibution  estimates  from  this  procedure  will  be  among 
the  inputs  to  the  retail  location  sector  of  the  Project  [7]. 
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II.  SOLUTION  METHODS  CONSIDERED 


D'Esopo  and  Lefkowitz  f 5]  describe  an  iterative  procedure  leading 
to  estimates  of  the  and  Oj  of  equation  (1).  This  algorithm  is  not  guaran¬ 
teed  to  always  converge,  and  hence  was  dismissed  as  unreliable. 

Wilson's  f 6 1  method  uses  one  of  the  best-known  of  the  "entropy- 
maximizimr"  approaches  to  travel  pattern  estimation.  The  two  principal 
objections  to  this  approach  were  (a)  its  extravagant  data  requirements, 
which  called  for  information  beyond  that  available  irom  the  census,  and 
(b)  the  unusual  and  nontraditional  form  of  the  function  g  of  equation  (1)  that 
results  from  the  algor ith  . 

Charnes,  Raike,  and  Bettinger  showed  (a)  that  the  t^  satisfying 
(1)  may  be  obtained  by  an  extremal  principle  of  convex  programming  form, 
and  (b)  that  this  principle  is  equivalent  to  a  standard  statistical  estimation 
procedure  of  the  information-theoretic  variety  [10]  --  a  hitherto  unknown 
fact.  The  CRB  procedure  has  none  of  the  lirr. nations  of  Wilson's  or  the 
D'Esopo-Lefkowitz  methods,  and  is  computationally  practicable.* 

III.  THE  SOLUTION  PROCEDURE 

Using  the  Kuhn-^ucker  optimality  conditions  Charnes,  Raike, 
and  Bettinger  [4]  pro\.-  the  following  Theorem.  Consider  the  convex 
programming  problem;  minimize  (t.^  In  K^t^  -  t^)  (2) 

^For  a  more  detailed  comparison  of  the  latter  two  solution  methods 
see  [8]  and  (12). 
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subject  to 


!  ■■■ '  D> 


i  ”1|  •  •  •  t  n 


S  t.  =  0{  j  =  l,  .  . .  ,  n 

3=1  lJ  1 

and  all  ty  *  0. 


If  a  travel  distribution  of  positive  t.  satisfies  the  gravity  equation  (1)  for 

some  set  of  potentials  t-  and  a ■  and  for  some  distance  function  K.  and  if 

1  J  >J 

these  ty  satisfy  the  origin-destination  requirements  (3),  then  these  t 
solve  the  convex  programming  problem  (2),  (3)  optimally.  Conversely, 
if  a  distribution  {t:  yl  solves  (2),  (3)  optimally,  then  potentials  t.,  cr 
exist  such  that  equation  (1)  is  satisfied  with  g(K^)  =  K.  . 

This  optimization  is  shown  to  yield  the  trip  distribution  which  is 
statistically  "least  distinguishable",  in  terms  of  the  statistical  information- 

9 

theoretic  measure  [4,  10] ,  from  the  ease-of-travel  distribution.  An 
additional  property  of  an  exact  optimal  solution  to  such  a  convex  program¬ 
ming  problem  ip  that  a  large  fraction  of  the  optimal  ty  must  be  non-integer. 
This  fact  has  gone  unmentioned  in  all  past  literature. 

Since  the  function  (2)  is  convex  and  separable  in  the  t.^,  it  may  be 
approximated  by  a  piecewise  linear  function  of  the  t  [1,  p.  348).  Thus 
transformed,  the  problem  would  take  the  form  of  a  "capacitated  distribution 
problem.  "  This  special  structure  can  be  solved  with  suitable  special 


Ease-of-travel  is  defined  as  the  inverse  of  the  generalized  travel 
cost,  normalized  to  a  0-1  interval. 
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computer  codes  two  orders  of  magnitude  faster  than  with  general  purpose 
linear  programming  codes.  Facing  squarely  the  fact  that  we  require 
integer  values  for  the  estimated  ty,  we  have  developed  a  new  two-piece 
approximation  which  guarantees  (for  extreme  point  solutions)  that  the  t^ 
will  be  integers. 

One  difficulty  in  the  capacitated  distribution  approach  lay  in  the 
fact  that  most  computer  codes  suitable  for  such  problems  will  accept 
only  integer  data.  Due  to  the  logarithmic  terms  in  the  functional  (2). 
truncating  the  coefficients  at  the  decimal  point  would  have  resulted  in 
catastrophic  rounding  errors.  Thus  it  was  necessary  to  scale  all  of 
the  input  data  upward  by  a  factor  of  10^  and,  of  course,  scale  the  resulting 
optimal  tjj  downward  by  the  same  factor. 

Initially  an  "out-of-kilter"  algorithm  code  was  employed  for 
solving  the  distribution  problem.  A  minor  complication  at  this  stage 
was  that  the  input  and  output  formats  of  the  out-of-kilter  code  were  not 
equipped  to  handle  integers  of  the  magnitude  of  the  "scaled-up"  input 
data.  Hence  these  formats  had  to  be  modified. 

The  all-in-core  program  was  run  in  422  seconds,  and  required  all 
of  the  core  memory  space  available  for  user  programs  on  the  CDC  6600 
at  the  University  of  Texas.  The  explanation  for  these  unusually  large 
time  and  core  requirements  may  lie  in  the  fact  that  (a)  all  districts  are 
available  as  destinations  for  each  district  of  origin,  and  (b)  the  2 -piece 
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linear  functional  necessitates  two  links  between  each  node  (district)  of 

the  network.  Thus  the  network  is  extremely  connected,  as  well  as  dense, 

compared  to  that  of  the  typical  capacitated  distribution  problem. 

Because  of  the  extraordinary  time  required  by  the  out-of-kilter 

algorithm,  the  problem  was  run  on  a  "primal"  algorithm  code  of  Dr. 

3 

Darwin  Klingman's.  An  optimal  solution  was  reached  in  only  nine  seconds 

A 

with  this  code.  Moreover,  the  core  requirement  was  reduced  by  three-fourths. 


IV.  DERIVATION  OF  THE  DIOPHANTINE  PIECEWISE  LINEAR  APPROXIMATION 
Recalling  the  nonlinear  equation  system  (2),  (3),  we  transform  the 


variables  t.^  as  follows: 


t. .  -*  ^ 


ij  iJ  U 


k=l 


(4) 


The  origin-destination  equations  (3)  become 


with  all  fi, .  >0. 
U  ~ 


n 

E 

ijl 

and 

n 

T 

3=1 

i*J 


£  6.  .k  =  D. 
k  =  l  IJ  3 


£  &  *  -  oi 

k=l  D  1 


j=l,  .  . . ,  n  (3a) 


i=l,  .  .  . ,  n  (3a) 


3 

We  wish  to  thank  Richard  Barr  of  the  Center  for  Cybernetic  Studies 
for  making  the  out-of-kilter  code  available  and  for  patiently  performing  the 
necessary  alterations,  and  Dr.  Darwin  Klingman  for  executing  the  primal 
code  calculation. 

4 Problems  much  larger  than  the  current  one  may  be  solved  by  any  of 
the  many  available  distribution  codes  which  are  not  executed  wholly  in-core. 
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The  functional  (2)  to  be  minimized  is 

s  (ttJ  m  kij  ty-tjj)  -$  (tij  In  “ij  +  ‘ij  ln  W 


(5) 


The  only  nonlinear  term  on  the  right-hand  side  of  (5)  is  of  form 
(t  ln  t).  Graphing  the  function  y=t  ln  t  (fig.  »: 


(The  function  has  a 


(  (e_1,  e"1) 

singularity  at  t=0;  the  derivative  is  infinite.  ) 

Figure  1 


It  is  this  function  that  must  he  approximated  by  a  piecewise  linear 
function.  Piecewise  linearity  is  achieved  by  introducing  upper  bounds 
on  the  6,/.  The  optima'  (extreme  point)  t(j  will  then  be  integral  or  non¬ 
integral  respectively  as  the  bounds  are  integral  or  non-integral.  The 
procedure  naturally  suggested  by  the  shape  of  the  cmve.  i.e..  to  place 
the  "knee"  of  the  approximation  at  (e'1.  e  ')  and  to  bound  the  6()  at  e 
is  therefore  useless  for  the  problem  at  hand,  because  the  resulting  optimal 

ty  will  not  be  integers. 

Bearing  in  mind  that  only  integer  results  are  competitive,  we  can 
mod  ty  the  functional  over  the  domain  of  its  non-integral  values.  We  intro¬ 
duce  the  modification  shown  in  Figure  2.  It  guarantees  integrality  and  an 
appropriate  functional  behavior  over  integer  values. 


0 


1 


t 


Figure  2 


Here,  the  slope  of  the  linear  function  is  -e  for  0  st  <1,  and  R  for 
1  <  t,  where  0  <  e  <<  1  and  R  =  min  0{,  max  D  }. 

1  J  J  ,  , 

1  1  <4 

Thus,  if  we  observe  the  condition  0*6*1,  then  -e6  +  R6  ~  t  In  t, 
where  t  =  6*  +  6^.  Therefore  the  functional  (5)  equals,  under  the  trans¬ 
formation  (4), 

-  •„*>  I"  K(J  ♦  (»„>  -  6ld2)  in  (6,/  ♦  #ij2)  -  («,/  ♦  •„*>! 


(6) 


«  +6i/'  *"  Kij  +  (-cSl/  +  RSU2)  -  (6i/  +  Sij2)*’ 

|j> 

which  is  clearly  linear  in  the  fiy  .  It  remains  to  minimize  the  right-hand 
side  of  equation  (6)  subject  to  the  conditions  (3a)  and  the  further  condition 
0  *  ^  *  1  for  all  l,  j. 


V.  RESULTS  AND  INTERPRETATIONS 

It  will  be  recalled  that  the  first  run  was  for  the  trip  type  "total 
internal  person  trips,  all  modes  of  travel"  (Figure  3).  We  anticipate 
further  applications  of  the  method  developed,  including  (a)  other  trip 
types  of  less  aggregated  character,  and  (b)  inclusion  of  the  major  traffic 
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arteries  ("external  stations")  entering  and  leaving  the  continuously 
urbanized  area. 

The  travel  distribution  estimate  obtained  from  the  procedure 
outlined  above  had  some  dominating  characteristics.  Traffic  volume 
was  concentrated  on  relatively  few  interdistrict  links.  Other  links  were 
assigned  token  transfers  many  orders  of  magnitude  less  than  those  of 
the  major  links.  Further,  a  pair  of  adjacent  districts  are  often  seen  to 
be  connected  by  a  major  link,  exchanging  approximately  equal  numbers 
of  travelers. 

A  potentially  important  result  of  this  study  is  the  observation 
that  the  distance  between  districts  (as  we  have  measured  it)  apparently 
exerted  strong  influence  on  the  distribution  estimate,  despite  the  fact 
that  in  equation  (1),  g(K.j)  =  Recall  that  in  many  previous  studies, 

an  effort  is  made  to  represent  g(K^)  as  K^a,  where  a  >  1,  giving  distance 
a  heavier  weight  in  the  estimate.  This  may  prove  to  have  been  a  super¬ 
fluous  endeavor. 

VI.  CONCLUSION 

Further  extension  of  the  work  of  Charnes,  Raike,  and  Bettinger  [4] 
is  in  progress,  including  a  new  direct  characterization  of  the  potentials 
Tj,  Oj  through  optimal  solution  of  an  unconstrained  convex  programming 
problem.  Corresponding  development  of  a  new  algorithm  for  the  t- 


follows.  These  characterizations  stem  from  the  dual  cf  a  variant  of  an 
extended  geometric  programming  problem  put  forward  in  another  context 
by  Charnes  and  Cooper  {2]. 
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